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Abstract 


We  report  on  the  verification  of  the  multiple  scattering  and  thermal  emission  model  PILOT-EX 
described  in  a  previous  publication  and  extend  the  scope  to  include  non-isothermal  clouds.  We 
also  verify  the  expressions  leading  to  the  thermal  version  of  the  radiative  transfer  equation,  fill  in 
the  details  of  the  Gaussian  spherical  cloud  formulation,  describe  methods  for  treating  the 
emissivity  correction  for  limited  bandpass  instruments,  and  discuss  the  significance  of  the  results 
in  terms  of  radiation  contrast  and  aerosol-induced  noise.  We  extend  the  model  to  include 
arbitrary  “embedded  sources”  and  verify  results  against  standards  in  the  literature  over  a  variety 
of  incident  angles  and  demonstrate  the  significance  of  the  results  in  computing  infrared  contrast 
for  use  in  imaging  systems  performance  analysis.  A  major  finding  is  that  comparisons  with 
known  (isothermal)  solutions  show  agreement  to  within  five  significant  figures  for  the  case  of 
plane  layers  for  both  normal  incidence  and  various  slant  path  angles.  Also,  preliminary  results 
for  the  case  of  temperature  stratification  demonstrate  a  significant  build  up  of  radiance  in  the 
(hotter)  cloud  interior,  which  leads  to  lower  values  of  the  apparent  emissivity  for  non-isothermal 
clouds. 
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Preface 


The  main  purpose  of  this  report  is  to  serve  as  a  follow  up  to  a  previous  report  describing  the 
main  features  and  utility  of  the  PILOT-EX  model  used  to  generate  exact  solutions  of  the  thermal 
and  optical  versions  of  the  radiative  transfer  equation.  There  are  roughly  seven  topics  covered, 
some  of  which  are  ancillary  to  our  previous  report  and  others  that  can  be  treated  as  “stand  alone” 
but  are  nevertheless  related  to  the  previous  work. 
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Executive  Summary 


This  report  is  the  second  in  a  series  of  reports  documenting  the  development  and  verification  of 
the  U.S.  Army  Research  Laboratory  PILOT-EX  aerosol  emissive  effects  model.  An  earlier  report 
covered  the  general  concepts,  a  short  theoretical  development,  and  some  preliminary  verification 
with  limiting  cases.  In  this  we  report  we  fill  in  many  of  the  mathematical  details  of  the  model, 
provide  many  more  cases  for  verification,  and  extend  the  scope  to  include  non-isothermal  clouds. 
We  also  verify  the  expressions  leading  to  the  thermal  version  of  the  radiative  transfer  equation, 
fill  in  the  details  of  the  Gaussian  cloud  formulation,  describe  methods  for  treating  the  emissivity 
correction  for  limited  bandpass  instruments,  and  discuss  the  significance  of  the  results  in  terms 
of  radiation  contrast  and  aerosol-induced  noise,  which  we  relate  to  the  problem  of  computing 
infrared  contrast,  familiar  to  the  user  community  and  used  in  imaging  systems  performance 
analysis.  On  the  theoretical  development,  we  also  extend  the  model  to  include  arbitrary 
“embedded”  radiation  sources  and  verify  results  against  standards  in  the  literature  over  a  variety 
of  incident  angles.  This  development  could  be  of  significance  to  modeling  (hot)  missile  plume 
signatures.  A  major  theoretical  finding  is  that  comparisons  with  known  (isothermal)  solutions 
show  agreement  to  within  five  significant  figures  for  the  case  of  plane  layers  for  both  normal 
incidence  and  various  slant  path  angles.  Also,  preliminary  results  for  the  case  of  temperature 
stratification  demonstrate  a  significant  build  up  of  radiance  in  the  (hotter)  cloud  interior,  which 
leads  to  lower  values  of  the  apparent  emissivity  for  non-isothermal  clouds.  The  details  given  here 
along  with  those  also  given  in  our  initial  report  offer  a  complete  solution  to  the  radiative  transfer 
equation  for  finite  aerosol  clouds  of  radial  symmetry.  The  method  addresses  both  aerosol  thermal 
emission  as  well  as  reflection,  including  all  orders  of  multiple  scattering.  The  method  represents 
exact  solutions  that  compare  favorably  with  other  studies  when  applied  to  plane  layers.  We  have 
also  demonstrated  the  significance  of  the  results  in  systems  evaluation  through  the  effect  on  the 
thermal  derivative,  or  “radiation  contrast,”  used  in  evaluating  infrared  scanning  systems  such  as 
the  FLIR  (forward-looking  infrared).  The  methodology,  per  se,  is  applicable  only  to  isothermal 
clouds  and  isotropic  scattering  although  these  are  not  necessarily  fundamental  limitations  on  the 
method.  In  yet  another  study  we  have  extended  the  method  to  include  non-isothermal  clouds, 
non-isotropic  scattering,  and  other  practical  effects  using  various  approximations. 


1.  Introduction 


The  main  purpose  of  this  report  is  to  serve  as  a  follow  up  to  a  previous  report,  which  described 
the  main  features  and  utility  of  the  PILOT-EX  model  used  to  generate  exact  solutions  of  the 
thermal  and  optical  versions  of  the  radiative  transfer  equation  [1].  There  are  roughly  seven  topics 
covered,  some  of  which  are  ancillary  to  our  previous  report  and  others  that  can  be  treated  as 
“stand  alone”  but  are  nevertheless  related  to  the  previous  work.  Sections  2, 4  and  5  are  examples 
of  the  former  type  which  cover  such  details  as;  the  formulation  of  the  radiative  transfer  equation 
to  include  natural  thermal  emission  in  section  2,  the  details  of  the  Gaussian  cloud  formulation 
borrowed  from  the  COMBIC  [2]  model  in  section  3,  and  the  verification  of  the  model  by 
comparisons  with  other  methods  from  the  literature  in  section  5.  In  section  4  we  give  an  example 
of  some  interesting  manifestations  of  the  exact  solutions  which,  when  applied  to  Gaussian 
clouds,  give  rise  to  various  effects  such  as  infrared  “edge  brightening,”  which  is  often  associated 
with  visible  wavelengths  and  natural  clouds.  We  extend  the  scope  of  the  study  somewhat  in 
sections  6,  7,  and  8,  by  addressing  such  issues  as  sensor  bandpass  effects  and  the  definition  of 
apparent  temperature  in  an  obscuring  environment  (section  6),  the  concept  of  radiation  contrast 
for  small  differential  signal  analysis  (section  7),  and  the  effect  of  aerosol-induced  path  radiance 
on  system  performance  parameters  such  as  the  signal-to-noise  ratio  (section  8).  In  section  9  we 
extend  the  formulation  to  include  thermal  stratification  and  define  a  practical,  or  “effective,” 
emissivity  to  augment  the  previous  work  which  treated  the  isothermal  case  only. 


2.  Thermal  Version  of  the  Radiative  Transfer  Equation 


Although  the  thermal  emission  term  discussed  in  section  3  of  our  previous  paper  [1]  seems  clear 
and  intuitive,  there  are  some  fundamental  difficulties  in  justifying  the  form  when  used  in  realistic 
atmospheric  applications.  The  root  of  the  difficulty  lies  in  the  use  of  the  Blackbody  function, 
B(X.,T),  which  requires  that  the  system,  atmosphere,  and  aerosol  cloud  be  in  strict 
thermodynamic  equilibrium,  a  condition  certainly  not  valid  for  any  realistic  scenario  that  we 
wish  to  address.  To  solve  the  non-equilibrium  problem  rigorously  would  require  an  in-depth  look 
at  the  kinetics  of  thermal  motion  and  going  far  beyond  the  relatively  simple  Maxwell-Boltzman 
treatment  in  the  classical  texts  and  certainly  beyond  the  scope  of  the  treatment  here.  Fortunately, 
there  is  a  practical  alternative  using  the  concept  of  “local  thermodynamic  equilibrium”  which  has 
been  rigorously  established  as  valid  for  most  atmospheric  applications  and  which  we  will  apply 
here  [3]. 

We  begin  with  the  differential  version  of  the  “standard”  radiative  transfer  equation  [4,5]  for 
computing  the  total  radiance,  I(r),  including  both  the  directly  transmitted  component  and  any 
contributions  to  multiple  in-scattering  as  (fig.  1): 

_L^=_(kabs  +osct)I(r)+osct  Jsrc(r)  (1) 
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where  C(r)  is  the  aerosol  mass  concentration  and  (kabs  >cSct)  are,  respectively,  the  mass  absorption 
and  mass  scattering  coefficients  which,  in  turn,  are  related  to  the  more  usual  mass  extinction 
coefficient  and  single  scattering  albedo  as 


<*  =  kabs+asc( 


(2a) 


“-=k 


abs  ^sca 


(2b) 


which  are  the  same  definitions  used  in  our  previous  paper,  although  we  have  altered  the  notation 
slightly  to  conform  with  cited  references.  The  major  unknown  in  eq.  (1)  is  the  so-named  (optical) 
source  function,  Jsrc(r),  which  we  ultimately  seek  to  determine.  The  differential,  dl(r),  in  eq.  (1) 
applies  to  a  small  cylindrical  volume  of  differential  thickness,  dr,  aligned  along  the  propagation 
vector,  r,  as  sketched  in  figure  (1). 


in-  scatter  I(r+dr) 


Figure  1.  Sketch  demonstrating  the  differential  radiative  transfer  equation. 

In  eq.  (1),  the  first  term  on  the  right  accounts  for  losses  due  to  both  absorption  and  out-scatter 
(i.e.,  extinction)  along  the  main  beam  and  the  second  accounts  for  in-scatter  of  ambient  radiative 
sources  (sky  and  surface)  from  all  directions.  Without  committing  to  any  specific  mathematical 
form,  we  next  add  an  ad  hoc  emission  term  to  the  right  side  of  eq.  (1)  as  follows: 


_L_dI(r) 
C(r)  dr 


=  -(kabs  +asct)I(r)+asct  Jsrc(r)+ 


_L^- 

C(r)  dr 


(3) 


where  the  added  term,  at  this  point,  could  include  any  type  internal,  or  “embedded,”  radiation 
source,  either  artificial  or  natural.  For  the  special  case  treating  thermal  emission,  it  is  plausible  to 
seek  an  expression  defining  an  aerosol  mass  emission  coefficient  analogous  to  the  mass 
absorption  and  mass  scattering  coefficients  such  that 


C(r)  dr 


=  eemsB(X,T) 


(4) 


where  B(X,T)  is  the  Planck  Blackbody  function,  X  is  the  wavelength,  and  T  is  the  real  (i.e., 
thermodynamic)  temperature  of  the  aerosol  and  can  be  a  function  of  r  if  we  include  thermal 
stratification  as  we  will  do  in  section  8.  In  eq.  (4)  we  have  introduced  eems,  as  an  ad  hoc  “mass 
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emission  coefficient”  analogous  to  the  mass  absorption  and  mass  scattering  coefficients,  kabS  and 
crsct.  To  specifically  address  the  case  of  natural  thermal  emission,  we  next  introduce  the 
equilibrium  version  of  the  radiative  transfer  equation  as  follows: 


1  dl 
C(r)  dr 


-kabsB(X,T)  + 


C(r)  dr 


(5) 


from  which  we  can  immediately  reproduce  a  version  of  Kirchoff  s  relationship  by  invoking  the 
equilibrium  condition  as  [dl/dr=0]  to  yield 


J_dU 

C(r)  dr 


=  kabsB(X,T) 


(6) 


which,  from  comparison  with  eq.  (4),  immediately  requires  the  mass  emission  coefficient  to  be 
equal  to  the  mass  absorption  coefficient  (i.e.,  eems=kabs),  which  is  the  essence  of  Kirchoff  s  Law 
as  applied  to  an  aerosol  medium. 

We  now  formally  apply  the  concept  of  local  thermal  equilibrium  to  the  problem  at  hand  by 
inserting  eq.  (6)  into  eq.  (3)  to  yield  the  following: 

=-(k*  +ara)I(r)+a„,  JUr)+k,l5B(X,n  P) 


which,  at  this  point,  we  cautiously  refer  to  as  an  approximation  and  make  explicit  the  fact  that 
the  (new)  source  function,  J'src,  must  now  include  not  only  the  effect  of  multiple  in-scatter  from 
external  sources  (as  before),  but  also  must  include  the  multiple  scattering  contribution  of  the 
internal  thermal  emissions  which  come  from  all  parts  of  the  cloud.  Next  we  rewrite  eq.  (7)  to 
explicitly  account  for  both  an  ambient  (external  in-scatter)  part  and  an  emissive  (internal 
scattering)  part  as 


1  dl(r) 
C(r)  dr 


=  ~(kabs  +  asct)I(r)+  asct  Jext(r)+<Wmt(r)  +  kabs  B(^T) 


(8) 


which  is  consistent  with  the  concept  defining  the  total  source  function  as  the  sum  of  the  ambient 
and  emissive  components  [i.e.,  J'src(  r)=  Jext(  r)  +  Jjnt(  r)]  as  was  done  our  earlier  work. 

The  degree  to  which  the  approximation  implicit  in  eq.  (8)  holds  depends  upon  the  degree  to 
which  the  atmosphere/aerosol  system  is  in  local  thermodynamic  equilibrium.  For  the  atmosphere 
in  general,  many  studies  have  shown  the  approximation  to  be  valid  up  to  an  altitude  of  about 
50,000  km,  above  which  the  Maxwell-Boltzman  distribution  is  not  necessarily  maintained  and, 
thus,  strict  equilibrium  as  applied  here  becomes  invalid.  The  inclusion  of  aerosol  particles  does 
not  present  a  big  conceptual  problem  because  it  is  the  surrounding  and  intervening  air  molecules 
that  actually  provide  the  thermal  agitation  required  to  maintain  thermal  equilibrium.  However, 
there  is  a  need  to  assure  good  thermal  mixing  so  that  the  temperature  of  the  aerosol  particles  and 


the  surrounding  air  molecules  remain  in  near  thermal  equilibrium,  which  is  a  valid  assumption 
here  if  we  assume  the  aerosol  particles  to  be  small.  There  may  be  other  applications  wherein  the 
accuracy  of  the  above  derivation  may  come  into  question;  however,  the  final  expressions  as 
given  here  have  been  more  or  less  been  accepted  as  valid  approximations  for  most  engineering 
applications  [6]. 

Owing  to  the  linearity  of  the  underlying  expressions,  it  is  clear  from  both  a  mathematical  and 
physical  sense  that  the  problem  can  be  solved  in  two  parts  by  separating  the  external  and  internal 
contributions.  That  is,  we  can  rewrite  eq.  (8)  in  the  form  of  two  differential  equations  that  we 
will  refer  to  as  the  “optical”  and  “thermal”  versions  of  the  radiative  transfer  equation.  That  is, 

optical  version : 

^  =  -IW  +  C0„J„,(T) 
dx 

thermal  version : 

-7-^  =  0  -  (t;T)  +  o),J  CO) 


where  we  have  also  made  two  major  substitutions;  the  first  by  using  eq.  (2)  to  insert  the  mass 
extinction  coefficient  [a=kabs+asca]  and  albedo  [ff>0=Gsca/(kabS+asca)],  and  the  second  to  replace  the 
spatial  variable,  r,  with  the  path  optical  thickness,  x  [=x(r,r0)].  Following  conventional  practice, 
we  now  define  the  optical  thickness  between  any  two  points  in  terms  of  the  concentration  path 
integral  as  [note  dx=aC(r)dr]: 


(9a) 

(9b) 


t(r„,r)=Jdx'=jaC(r-)dr’  (10) 

where  it  is  understood  that  the  integration  proceeds  along  the  straight  line  path  connecting  the 
two  points  r0  and  r. 

For  the  special  case,  Jext  =Jint  =0,  both  expressions  in  eq.  (9)  can  be  integrated  immediately  to 
yield  the  following  “zeroth”  order,  or  Beer’s  Law,  solutions;  that  is, 

optical  version : 
thermal  version : 

Jdir(x)  =  (l-a)0)jB*AX(x',T)e-^'r»)dx' 


(11a) 

(Hb) 


which  are  identical  in  form  to  what  we  termed  the  direct  contributions  discussed  in  [1]. 
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The  differential  expressions  of  eq.  9  can  be  converted  to  integral  expressions  by  formally 
integrating  both  sides  of  the  equations.  Following  the  usual  procedure,  we  first  multiply  each 
side  by  an  integrating  factor  of  the  form  e4^  and  then  integrate  along  an  arbitrary  path  whereupon 
we  have,  from  eq.  (9a),  after  some  slight  rearrangement, 


+  I(i  V  Uf  =  <o„  Jj„  (i')e'di' 


0 

or 


(12a) 


(12b) 


where,  in  the  second  expression,  the  term  on  the  left  side  follows  directly  from  the  chain  rule  for 
differentiation  and  can  be  integrated  immediately  to  yield 


optical  version : 

I(i)  =  I0e-’+BIJ'j„,(i')e-|’-’',di' 


(13a) 


where  we  have  performed  the  (trivial)  integration  on  the  left  side  of  eq.  (12b)  and  divided 
through  by  the  integrating  factor,  eT .  Written  as  such,  eq.  (13a)  is  the  integral  form  equivalent  of 
eq.  (9a).  The  treatment  for  the  thermal  version  proceeds  in  the  same  manner  and  ultimately  leads 
to 


thermal  version : 

I(T)  =  0  -  )  '}b(X,  T)e"<t-T  ’dr'  +  <o„  Jjin,  (i')e~<T-1  Mi' 


(13b) 


which  accounts  for  all  radiance  emanating  from  the  aerosol  cloud  and  does  not  include  the  direct 
contribution  from  the  target  which  has  already  been  accounted  for  in  the  “optical”  version.  We 
can  think  of  the  first  term  on  the  right  side  of  eq.  (13b)  as  a  direct  (emissive)  component 
accounting  for  all  direct  (aerosol)  emissions  along  the  path,  very  much  analogous  to  the  (hard 
target)  direct  term  of  eq.  (13a).  Equations  (13  a)  and  (13b)  are  valid  for  any  geometric  shape  and 
form  the  basis  for  the  multiple  scattering  solutions  throughout  our  previous  work  in  [1]. 

Although  we  have  focused  here  on  (natural)  thermal  emission,  it  turns  out  that  the  derivation  can 
also  be  applied  to  the  more  general  case  treating  any  (artificial)  embedded  source  of  arbitrary 
form.  The  derivation  is  similar  to  that  described  above  except  that  we  need  not  necessarily  be 
concerned  about  the  condition  of  thermal  equilibrium.  Thus  we  recast  eq.  (4)  in  a  more  general 
form  as 


1  dl 


emb  _ 


C(r)  dr 


JembOO 


(14) 
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where  jemb  is  referred  to  as  the  “emission  coefficient”  and  may  or  may  not  be  a  function  of  the 
cloud  temperature.  In  eq.  (14)  we  have  used  the  subscript  “emb”  to  stand  for  “embedded” 
sources  following  the  general  ideas  of  the  classical  texts.  We  can  now  retrace  the  steps  leading  to 
eq.  (13)  and  introduce  one  more  expression  to  address  the  general  emissive  case  as 

embedded  sources : 

ICO  =  J-Lb('Oe~(T‘t')d'c'  +  co0  {jemb(t')e"(t_t,)dT'  ^ 


where  jemb(x)  is  some,  as  yet  unspecified,  function  serving  the  same  purpose  as  the  direct 
emissive  component  in  the  thermal  case  and  Jemb(0  is  the  multiple  scattering  contribution.  In 
section  4  we  will  test  eq.  (15)  and  compare  with  the  classic  results  of  Van  de  Hulst  [5]. 


3.  Gaussian  Cloud  Formulation 


In  this  section  we  develop  the  “standard”  Gaussian  cloud  formulation  used  in  the  source  function 
calculations.  The  main  objective  is  to  develop  the  analytical  expressions  needed  for  calculating 
the  optical  thickness,  x(r  r0) ,  between  any  two  points,  rand  r0,  assuming  a  Gaussian  aerosol 
cloud  centered  at  the  origin  as  sketched  in  figure  2.  With  reference  to  figure  2,  the  task  is  to 
calculate  the  optical  thickness  along  the  dashed  path  beginning  at  some  arbitrary  point  r  and 
proceeding  outward  to  the  practical  “edge”  of  the  cloud  indicated  by  the  dashed  arc  at  r0,  which 
corresponds  to  Ro,  the  cloud  “practical”  radius. 


Figure  2.  Sketch  of  Gaussian  cloud  formulation  and  nomenclature. 

For  a  constant  mass  extinction  coefficient,  a,  it  follows  from  the  definition  of  eq.  (10)  in  the 
previous  section  that  the  major  computational  burden  lies  in  computing  the  path  integrated 
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concentration,  or  “CL  Product,”  which,  for  the  case  of  a  Gaussian  concentration  distribution,  is 
written  as 


CL(r,r  )  =  jc(r')dr'  =-£'■-  Jexp-f-^dr' 
,  (V2tig0)3  1  V2a0 


(16) 


where,  as  noted  in  figure  2,  the  point  indicated  by  the  vector,  r'  also  corresponds  to  a  vector 
distance  x'  from  the  integration  starting  point. 

From  simple  vector  geometry  and  figure  2  we  have  [r'=r+x'],  which  we  next  substitute  into 
eq.  (16)  to  yield 


CL(r,r0) 


(17) 


where  we  have  also  made  explicit  the  fact  that  the  integration  proceeds  along  the  straight  line 
path,  x',  starting  at  x'=0  and  ending  at  x'=x0,  as  indicated  in  the  sketch  of  figure  2.  In  eq.  (17), 
the  parameters  Q0  and  a0  are  known  constants  representing,  respectively,  the  cloud  total  mass 
(grams)  and  the  Gaussian  half  width  (meters)  and  are  determined  in  a  manner  described  later  in 
this  section. 

The  next  objective  is  to  express  eq.  (17)  into  a  more  workable  scalar  form.  We  do  this  by  first 
expanding  the  exponential  argument  in  the  above  expressions  by  applying  simple  vector 
geometry  to  obtain  the  following: 

r'2  =[r +  x']2 

=|r|2  +|x'|2  +2rx'(f -x') 

=  x'2 +r2 +2rx'p 
=  (x'  +  pr)2  +  r2(l-p2) 


where  p  is  the  cosine  of  the  angle  between  the  radial  vector,  r,  and  the  propagation  vector,  x'  (or 
x0 ,  since  the  two  are  collinear).  The  form  of  the  last  expression  in  eq.  (18)  follows  from  simple 
geometry  with  some  minor  manipulation  and  can  be  verified  by  direct  expansion. 

We  now  rewrite  eq.  (17)  with  the  aid  of  eq.  (18)  into  a  more  workable  (scalar),  form  for 
subsequent  integration  as  follows: 
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where  we  have  made  explicit  the  fact  that  the  expression  as  now  written  is  dependent  upon  the 
scalar  variables  (r,r0,p).  The  constant,  Q  ,  turns  out  to  be  numerically  equal  to  aerosol  mass 
concentration  at  the  integration  starting  point  [i.e.,  Ci=C(r)]  and  can  be  calculated  accordingly.  It 
is  of  some  importance  to  note  from  the  form  of  eq.  (19)  that  the  final  integrated  result  is 
dependent  only  upon  the  relative  angle,  p,  and  the  scalar  distances  (r,Xo)  where  x0  is  the  total  path 
distance  [x02=r2+r02].  This  fact  was  found  to  be  important  in  applying  the  mathematics  for 
evaluating  the  source  function  in  [1]. 

Following  the  usual  approach  we  next  introduce  the  intermediate  variable,  t'2=(x'+pr)2/2o02, 
allowing  us  to  compress  eq.  (19)  into  the  following  form: 

*o 

CL(r,,r2;p)  =  V2o0C,  Jexp[-t'2]dt'  (20) 


where  we  have  also  substituted,  dt'=dx72cr0.  The  major  advantage  gained  from  all  this 
manipulation  lies  in  the  fact  that  the  integral  in  eq.  (20)  can  now  be  expressed  in  terms  of  the 
well-known  error  function  and  can  thus  be  integrated  immediately  as  follows: 


CL(r,r0 ; p)  =  V2a0C,  { Jexp[  - t'2  ]dt'  -  JexpR'2  ]dt'} 

0  O 

=  Vrc/ic^C,  (erf(t0)  -  erf (t)} 


=  CL  0{erf(- 


-)-erf(r 


(21) 


where  CL0  is  a  another  constant  ( CL0  =  'Jn/2a0Cl )  which  turns  out  to  be  numerically  equal  to 
the  total  “center-to-edge,”  or  CL  radial  path  and  will  be  discussed  more  fully  in  the  following 
paragraphs. 

The  error  function  as  used  in  eq.  (21)  is  defined  formally  for  our  applications  in  the  usual 
notation  as 
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(22) 


which  can  be  evaluated  readily  using  any  of  a  number  of  numerical  approximations  described  in 
the  standard  handbooks  [7]. 

In  actually  applying  the  above  formulation  to  the  task  at  hand,  all  calculations  are  performed  on  a 
“standard”  cloud  of  fixed  radius  Ro  representing  the  “practical”  extent  of  the  Gaussian  cloud 
(which  strictly  extends  to  infinity).  For  practical  purposes  we  arbitrarily  define  the  edge  as  the 
point  where  the  concentration  drops  to  1/1000  of  the  value  at  the  centroid.  We  furthermore 
define  the  total  mass  of  this  standard  cloud  such  that  the  total  edge-to-edge  CL  path,  or  “CL 
diameter,”  calculates  to  unity.  After  some  manipulation,  this  leads  to  the  following  values  for  the 
“standard”  cloud: 


cr0  =  R0  /^/2Ln(1000)  =13.452 
Q0  =  2tcg2  /(erf  (R0  /V2g)  =  1350.0 


(23) 


where  the  numerical  values  shown  correspond  to  calculations  performed  using  a  cloud  of  one 
meter  diameter  [i.e.,  Ro=Do/2=0.50(m)]  and  are  the  actual  values  used  throughout  the 
calculations.  With  the  Gaussian  parameters  so  determined,  the  calculation  of  the  CL  Product 
over  any  path  inside  the  cloud  is  readily  found  from  eq.(21)  from  which  the  optical  thickness  is 
then  simply 


x(r,r0;p)  =  aCL(r,r0;p) 


(24) 


which  follows  directly  from  the  defining  equation,  eq.  (10),  of  the  previous  section. 

It  is  important  to  note  from  eq.  (24)  that  we  can  vary  the  cloud  optical  thickness  in  either  of  two 
ways:  (1)  either  by  changing  the  optical  variable,  a,  or  (2)  by  changing  the  cloud  concentration 
through  the  Gaussian  parameters  Q0  or  a0.  It  turns  out  that  these  two  methods  are  totally 
equivalent  and  are  part  of  the  basis  for  the  optical  thickness  (not  r)  being  the  “real”  independent 
variable  for  the  optical  calculations.  This  means  that  once  we  have  done  the  calculations  for  our 
“standard”  unit  cloud,  the  results  are  valid  for  any  cloud  of  the  same  optical  depth,  regardless  of 
the  spatial  extent — a  circumstance  that  makes  the  formulation  much  more  universal  than  it  might 
at  first  appear. 

To  summarize  this  section  we  show  a  numerical  example  in  figure  3  of  some  edge-to-edge  cross- 
sectional  plots  of  concentration  and  CL  Product  for  our  standard  cloud  using  the  parameters  of 
eq.  (23)  assuming  a  unit  cloud  diameter  (i.e.,  Ro=Do/2=0.50). 
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Figure  3.  Example  of  optical  and  spatial  representation  for  Gaussian  cloud. 


Referring  to  the  upper  two  plots  of  figure  3,  those  on  the  left  are  based  upon  the  usual  spatial 
representation  where  the  abscissa  represents  linear  distance  referenced  from  the  (left)  edge  of  the 
cloud  and  passing  and  through  the  cloud  centroid.  The  plots  on  the  right  represent  the  same 
numerical  data,  but  we  use  the  CL  Product  as  the  abscissa  instead  of  the  spatial  distance.  Note 
that  we  label  the  axes  as  either  normalized  spatial  distance  (d/do)  or  normalized  optical  distance 
(x/x0).  In  both  cases  the  curves  labeled  (a)  represent  the  concentration  (normalized  to  unity  at  the 
center),  the  curves  labeled  (b)  represent  the  CL  Product  referenced  to  zero  from  the  left  edge,  and 
the  curves  labeled  (c)  represent  the  “radial”  CL  Product  reference  to  zero  at  the  midpoint.  Note 
from  the  plots  the  obvious  radial  symmetry  with  the  concentration  maximum  at  the  center  and 
the  CL  Product  rising  from  a  value  of  zero  at  the  left  edge  to  a  value  of  0.50  at  the  center  and 
finally  to  a  value  of  unity  at  the  right  edge.  The  radial  CL  Product  is  zero  at  the  center  and 
increases  to  a  maximum  value  of  0.50  toward  either  edge.  Clearly,  in  the  plots  on  the  right  the 
optical  thickness  appears  as  a  straight  line  as  it  must  and  the  concentration  plot  is  much  broader 
in  appearance. 

In  the  lower  plots  of  figure  3,  we  show  another  example  in  how  the  choice  of  representation 
affects  the  appearance  of  the  source  functions  described  in  [1].  The  particular  examples  shown 
here  represent  the  emissive  source  function  based  on  a  cloud  optical  thickness,  or  “optical 
diameter,”  of  unity  (i.e.,  x0=l).  The  various  curves  in  each  plot  represent  albedo  ranging  from 
zero  and  increasing  to  unity  in  the  direction  shown.  Clearly,  the  main  difference  between  the  two 
representations  is  that  the  optical  version  produces  generally  broader  plots  than  the  spatial 
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versions  and  might  offer  some  analytical  advantages.  Our  experience  is  that  the  optical 
representation  is  more  natural  and  is  usually  easier  to  work  with  analytically. 


4.  Radiance  Distribution  in  Gaussian  Clouds 


In  our  main  paper  we  made  mention  that  the  behavior  of  the  reflectivity  functions  could  lead  to 
some  interesting  effects  on  cloud  images.  In  particular  the  fact  that  the  reflectivity  functions 
show  a  maximum  at  moderate  optical  thickness  can  lead  to  an  edge  brightening  under  certain 
conditions  that  is  often  observed  for  natural  clouds  at  visible  wavelengths.  Some  insight  into  how 
this  comes  about  can  be  gained  by  looking  more  closely  at  the  various  relationships  existing 
among  the  transmissivity,  taer,  emissivity,  8aer,  and  reflectivity,  raer,  functions,  especially  the 
continuity  relationship,  which  we  rewrite  here  as 

tarcCO  +  Saerto  +  raerC^1  (25) 


which  is  valid  for  any  path  through  an  isothermal  cloud. 

We  demonstrate  these  relationships  in  the  two-dimensional  arrays  of  figure  4  which  were 
generated  by  applying  eq.  (30)  of  our  main  paper  along  multiple  parallel  paths  traversing  the 
entirety  of  a  Gaussian  cloud.  In  constructing  figure  4,  we  assumed  a  100  by  100  m  Gaussian 
array  of  1  m  resolution  and  standard  deviation,  a0,  of  30.0  m.  In  figure  4  the  first  three  columns 
refer  to  transmissivity,  emissivity,  and  reflectivity,  respectively,  as  calculated  using  the  source 
functions  for  an  isothermal  cloud.  Each  row  in  figure  4  corresponds  to  a  specific  edge-to-edge 
centerline  optical  thickness  ranging  from  x=0.50  (upper  row)  to  r=8.0  (bottom  row)  as  indicated 
in  the  side  margins.  Out  of  curiosity  we  have  also  included,  in  the  last  column,  a  difference  array 
between  the  emissivity  and  reflectivity,  which  according  to  the  continuity  condition,  should  be  of 
the  form  [l-e'T],  In  all  cases  the  albedo  is  set  to  co0=0.50. 
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Figure  4.  Two  dimensional  arrays  of  transmissivity,  emissivity  and  reflectivity. 


In  analyzing  the  results  of  figure  4,  the  transmissivity  (first  column)  and  emissivty  (second 
column)  trends  are  the  most  easy  to  interpret  in  that  they  show  a  general  monotonic  increase  in 
magnitude  with  increasing  optical  thickness,  which  also  gives  rise  to  an  apparent  increase  in  the 
perceived  “size”  of  the  resultant  cloud.  The  reflectivity  (third  column)  shows  a  similar  increase 
with  increasing  optical  thickness  but  not  nearly  as  marked  as  is  the  case  for  emissivity.  This 
difference  in  cloud  “brightness”  occurs  in  spite  of  the  fact  that  the  (single  scattering)  albedo  of 
0.50  indicates  equal  fractions  of  absorption  and  scattering  at  the  particle  level.  The  most  curious 
feature  in  figure  4  is  the  appearance  of  a  darkened  center  for  the  reflective  case  that  becomes 
more  and  more  pronounced  as  the  optical  thickness  increases  and  is  perceived  as  an  actual 
decrease  in  reflectivity  at  the  center.  This  “darkening”  effect  is  real  and  is  related  to  the 
maximum  in  reflectivity  discussed  in  our  main  paper,  but  is  more  easily  explained  as  a 
manifestation  of  the  continuity  relationship. 


We  probe  this  “brightening”  effect  more  quantitatively  in  figure  5  where  we  show  side-by-side 
cross-sectional  plots  of  emissivity  and  reflectivity  for  a  wide  range  of  situations. 
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Figure  5.  Cross-sectional  plots  of  transmissvity,  emissivity,  and  reflectivity. 


The  plots  in  figure  5  show  results  of  calculations  for  both  emissivity  (left)  and  reflectivity  (right) 
for  three  values  of  albedo:  .90  (top),  0.50  (middle),  and  0.10  (bottom).  Each  subplot  includes  1 1 
curves  representing  optical  thicknesses  ranging  from  0.125  to  10.0  which  occur  in  an  increasing 
monotonic  order  with  the  exception  for  the  center  values  in  the  reflective  case.  In  figure  5.  the 
general  trend  for  the  emissivity  curves  is  toward  a  monotonic  increase  in  emissivity  with 
increasing  optical  thickness  which  holds  for  all  albedos.  It  is  also  generally  true  from  the  three 
examples  that  the  effect  of  an  increased  albedo  is  to  lower  the  emissivity  (note  y-axis  scales) . 
There  is  also  a  tendency  for  the  emissivity  curves  to  “flatten  out”  and  mesh  together  as  the 
optical  thickness  increases,  an  effect  that  is  more  and  more  pronounced  as  the  albedo  decreases, 
and  in  the  extreme  gives  rise  to  a  broad  maximum  in  the  curves  which  mesh  toward  unity  at  the 
highest  optical  thicknesses.  The  reflectivity  curves  show  a  similar  behavior  except  that  the  effect 
of  albedo  is  reversed  in  the  sense  that  an  increased  albedo  gives  rise  to  a  higher  reflectivity  (note 
y-axis  scale  changes),  and  the  curves  not  only  tend  to  “flatten  out”  but  actually  decrease  in  the 
center  as  the  optical  thickness  increases  beyond  a  value  between  about  1  or  2.  From  figure  5,  the 
effect  becomes  more  and  more  marked  as  the  albedo  decreases.  This  decrease  in  the  center,  also 
evident  in  figure  4,  is  often  perceived  as  an  increase  at  the  cloud  edges  often  observed  in  natural 
clouds  in  the  sunlit  sky  and  has  also  been  observed  in  single  scattering  models  at  visible 
wavelengths  [8]. 
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5.  Verification  of  Pilot-Ex  (Isothermal  Plane  Layers) 


Although  the  main  emphasis  here  is  on  finite  clouds  of  spherical  symmetry,  it  turns  out  that  the 
multiple  scattering  method  described  in  our  earlier  paper  also  applies  to  plane  layer  symmetry  for 
which  there  are  other  studies  to  compare  results,  the  most  comprehensive  being  the  extensive 
tables  published  by  van  de  Hulst  [5].  Our  purpose  in  this  section  is  to  use  the  plane  layer  results 
as  a  verification  of  our  method  for  computing  the  emissive  and  reflective  source  functions.  As  it 
happens,  there  is  actually  very  little  data  for  direct  comparisons  of  the  source  functions,  per  se; 
however,  there  are  extensive  tables  for  the  integral  quantities,  like  reflection  and  to  a  lesser 
extent  emission,  which  will  be  main  basis  for  the  comparisons.  A  sketch  of  the  relevant  geometry 
for  the  plane  layer  is  shown  in  figure  6. 


Figure  6.  Sketch  demonstration  Plane  Layer  nomenclature. 


For  the  plane  layer  case  we  need  only  replace  the  three-dimensional  Gaussian  concentration 
functions  with  the  more  easily  treated  one-dimensional  form,  which  ultimately  results  in  the 
following  simple  expression  for  the  slant  path  optical  thickness: 

T(z,n)  =  ^  (26) 

u 

where  tv(z)  [=t(z,1)]  is  the  optical  depth  measured  over  the  vertical  spatial  distance,  z,  and  u 
[u=|p|],  is  the  absolute  value  of  the  cosine  of  the  slant  path  angle  (fig.  6).  For  this  study  the 
vertical  optical  depth  is  modeled  as  a  Gaussian  function  with  the  centroid  at  the  vertical  midpoint 
of  the  layer  and  again  is  based  on  a  “standard”  cloud  with  parameters  selected  such  that  the 
vertical  edge-to-edge  CL  product  taken  through  the  entirety  of  the  layer  calculates  to  unity. 

It  follows  directly  from  inspection  of  eq.  (26)  that  the  condition  of  azimuthal  symmetry  holds  so 
that  we  can  apply  the  methodology  of  section  3  in  our  previous  report.  It  also  turns  out  that  for 
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the  plane  layer  case  the  (numerous)  angular  integrations  can  be  performed  analytically  by 
employing  the  Exponential  Integral, 


+1 

1=  Je't/udu  =  E2(T/u)  (27) 

where  Ea(x)  is  the  well  known  second  order  Exponential  Integral  which  was  evaluated  using 
established  numerical  approximations  [7].  Other  than  Eq.  (26)  and  (27)  there  were  no  other 
substantive  changes  needed  to  address  the  plane  layer  case. 

In  figure  7  we  show  the  results  of  calculations  for  the  plane  layer  using  the  same  general  scheme 
that  was  used  for  the  spherical  clouds  in  our  earlier  work  (fig.  5)  [1].  The  calculations  are  based 
upon  the  normalized  source  functions  and  thus  vary  in  overall  magnitude  from  zero  to  unity.  The 
abscissa,  for  the  plane  layer  case,  is  the  edge-to-edge,  or  better,  “top-to-bottom”  normalized 
optical  thickness  referenced  along  the  vertical  from  top  to  bottom. 


Figure  7.  Source  functions  for  plane  layers:  (a)  reflection,  (b)  emission. 

The  layout  in  figure  7  is  identical  to  that  of  figure  3  and  show  the  reflective  case  on  the  left  and 
the  emissive  case  on  the  right.  For  each  case  there  are  six  subplots  corresponding  to  six  values  of 
cloud  total  optical  thickness  ranging  from  xo=0.250  to  8.0,  as  indicated  on  the  figures.  The 
various  curves  in  each  subplot  correspond  to  albedos  ranging  from  zero  to  unity  in  steps  of 
Acd=0.10  and,  as  before,  either  increase  in  overall  magnitude  with  increasing  albedo  (reflective) 
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or  decrease  in  overall  magnitude  with  increasing  albedo  (emissive),  as  indicated  by  the  direction 
of  the  arrow  on  the  uppermost  plots.  The  calculations  were  performed  over  a  vertical  path  (i.e., 
p=l)  so  that  the  comparisons  with  the  spherical  results  (taken  along  the  cloud  centerline)  are  on  a 
common  basis  as  far  as  the  total  path  optical  thickness  is  the  concern. 

In  comparing  with  the  spherical  case,  we  note  from  figure  7,  the  same  general  “flattening”  of  the 
source  functions  for  lower  optical  thickness,  and  conversely,  the  increased  curvature  as  the 
optical  thickness  increases.  In  all  cases,  the  reflective  curves  are  concave  upward  and  the 
emissive  curves  are  concave  downward  as  they  also  are  for  the  spherical  case.  It  is  evident  from 
comparison  with  figure  5  and  our  earlier  work  in  [1]  that  there  are  no  qualitative  differences 
between  the  plane  layer  solutions  and  the  spherical  solutions.  In  fact,  there  are  very  little 
quantitative  differences  at  small  optical  thickness  where  all  solutions  converge  to  the  first  order 
result  (i.e.,  curves  are  flat).  At  the  higher  optical  thickness,  however,  the  plane  layer  curves  are 
significantly  lower  than  for  the  spherical  case  as  one  might  intuitively  expect. 

It  turns  out  that  there  are  relatively  few  cases  where  we  can  directly  compare  the  source 
functions,  per  se,  with  results  from  the  literature;  however,  there  are  extensive  tables  for 
comparing  the  integral  properties  such  as  reflectivity  which  we  show  in  table  1  for  the  special 
case  of  normal  incidence  (i.e.,  p=T).  The  comparisons  in  table  1  all  address  the  plane  layer  case 
using  results  from  either  the  methodology  of  our  previous  work  (labeled  Ref.  1)  or  results  from 
table  12  of  the  Van  de  Hulst  text  (labeled  Ref.  4). 

Table  1.  Comparison  of  calculated  reflectivities  for  normal  incidence  with  plane  layers. 


REFLECTION 


E35H 

wo=0.99 

wo=0.80 

wo=0.40 

wo=0.20 

Tau 

Ref.  1 

Ref.  5 

Ref.  1 

Ref.  5 

Ref.  1 

Ref.  5 

Ref.  1 

Ref.  5 

.03125 

.03044 

.03054 

.  02426 

.02426 

.  01178 

.01178 

.  00581 

.00582 

.06250 

.05990 

.05990 

. 04723 

.04723 

.02246 

.  02246 

.  01096 

. 01097 

. 12500 

.11606 

.  11605 

. 08978 

.  08978 

.  04119 

.  04119 

.  01978 

.01978 

.25000 

.21806 

.21806 

.  16305 

. 16305 

.07046 

.  07046 

.03299 

.03295 

.50000 

.38644 

.38581 

.27139 

.27139 

.10648 

.  10647 

.04808 

.04808 

1.00000 

.61602 

.61602 

.38439 

.38437 

.13056 

.13055 

.05637 

.05636 

2.00000 

.82786 

.82785 

.41663 

.41655 

.11971 

.11966 

. 04980 

. 04977 

4.00000 

.90228 

.90221 

.34105 

.34079 

.09234 

.09221 

.03847 

.03841 

8.00000 

.85092 

.85148 

.29033 

.28956 

.08403 

.08364 

.03551 

.03531 

Before  proceeding  there  are  some  potentially  confusing  circumstances  that  need  to  be  cleared  up 
when  comparing  with  the  Van  de  Hulst  results.  That  is,  in  the  Van  de  Hulst  text,  the  “reflection” 
function,  as  we  define  it  here,  is  separated  into  two  components:  one  labeled  transmission, 

T(t,u),  for  general  downward  (sky-to-surface)  diffuse  propagation  and  another  labeled  reflection, 
R(x,u),  for  general  upward  (surface-to-sky)  diffuse  propagation.  This  is  a  convention  dating  from 
early  astrophysical  studies  based  on  planetary  plane  layers  as  viewed  inward.  Furthermore,  the 
Van  de  Hulst  transmission  function  includes  a  direct  component,  e‘T ,  that  we  need  to  subtract  out 
to  make  the  comparisons  valid.  One  needs  also  to  be  cautioned  that  all  results  here  pertain  to 
Lambertian  sources  over  the  full  4n  upper  and  lower  hemispheres  as  the  driving  terms.  That  is, 
the  relationship  between  our  reflectivity,  raer,  and  the  Van  de  Hulst  functions  are  as  follows: 


raer(x,u)  =  R(t,  u)  +  T(t,u)  -  e'1 


(28) 
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where  the  R(x,p)  and  T(x,p  )  values  are  found  under  the  heading  of  “Lambertian  sources”  in 
table  12  of  the  Van  de  Hulst  text.  We  also  note  that  the  Van  de  Hulst  tables  refer  to  “intensities” 
referenced  outward  as  calculated  assuming  unit  input  and  are  equivalent  to  our  “normalized” 
source  functions.  The  results  here  all  pertain  only  to  reflection  since  this  is  the  only  case  reported 
in  the  Van  de  Hulst  tables;  however,  there  are  some  indirect  graphical  data  for  emission  which 
we  address  later  in  this  section. 

Looking  now  at  the  numerical  results  of  table  1,  we  see  clearly  that  the  comparisons  between  our 
work  and  Van  de  Hulst  are  quite  favorable,  showing  discrepancies  of  less  that  0.10  percent  in 
almost  all  cases  with  the  larger  exceptions  occurring  only  at  the  very  highest  optical  depths 
where  multiple  scattering  effects  are  most  significant.  Later  studies  verified  that  these  disparities 
could  be  further  diminished  by  using  a  higher  spatial  resolution  in  the  numerical  integrations, 
using  a  lower  threshold  for  convergence,  or  both.  Overall  the  results  here  are  qualitatively  similar 
to  the  spherical  case,  as  one  should  expect.  In  particular,  looking  carefully  at  the  numerical 
results  of  table  1  we  see  that  all  values  tend  to  monotonically  increase  and  reach  a  maximum 
somewhere  between  x0=2-4  beyond  which  the  values  again  decrease.  This  peak  in  the  reflection 
curve  is  very  similar  to  that  found  for  the  spherical  case  discussed  in  the  previous  section. 

In  passing  ,  it  also  needs  to  be  pointed  out  that  the  use  of  the  analytical  form  of  eq.  (27)  to 
perform  the  angular  integrals  was  more  of  a  necessity  than  a  simple  convenience  because  it 
turned  out  that  the  use  of  “brute  force”  numerical  methods  were  extremely  inefficient  and  in 
many  cases  required  angular  meshes  on  the  order  of  milliradians  to  obtain  the  highest  accuracies. 
This  problem  was  of  great  practical  concern  and  was  completely  resolved  by  using  the 
Exponential  Integrals  of  eq.  (27).  In  any  case,  based  upon  the  comparisons  here  there  seems  to  be 
little  doubt  as  to  the  validity  of  our  approach — at  least  on  matters  pertaining  to  reflection. 

In  addressing  the  emissive  case,  Van  de  Hulst  has  taken  the  related  approach,  which  we  referred 
to  as  “embedded  sources”  in  the  discussion  leading  to  eq.  (15)  in  section  1.  In  this  case  there  are 
no  extensive  tables  with  which  to  compare  our  results;  however,  Van  de  Hulst  does  give  one 
graphical  example  demonstrating  the  angular  dependence  which  we  can  relate  to  our  case  of 
reflection  and  emission.  The  relevant  expressions  for  comparison  in  the  plane  layer  case  are  as 
follows: 


Iref(X>U)  =  Pref[T05®o5T,(Z)]e'T(Z)/Udx' 

o 

Iemb(x>u)  =  Jjemb[T0,co0;T'(z)]e'T'(z)/udT' 


(29a) 

(29b) 


where  u  is  the  absolute  value  of  the  viewing  angle  cosine  and  z  is  vertical  distance  (fig.  6).  The 
source  functions  are  normalized  to  unity  for  the  case  of  reflectivity  and,  equivalently,  to  a 
uniform  unit  embedded  source  for  the  case  of  emission  as  described  by  Van  de  Hulst. 


17 


1  REFLECTION 

0.8 
0.6 
0.4 
0.2 
0 

Optical  Distance 

1 

0.8 
0.6 
0.4 
0.2 
0 

0  -  0,5  1 

Zenith  Cosine,  u 


0  0.5  1 

Optical  Distance 

EMISSION 


Zenith  Cosine,  u 


Figure  8.  Relative  intensity  functions  for  plane  layers:  (a)  reflectivity,  (b)  emission. 


Results  of  comparisons  for  a  single  case  (x0=T)  are  shown  in  figure  8.  The  upper  curves  in 
figure  8  refer  to  the  normalized  source  functions  and  are  added  here  only  for  convenience  to  the 
reader;  however,  it  should  be  noted  that  the  emission  functions  here  are  not  normalized  to  unity 
and  that  the  ordering  of  the  emissive  curves  are  not  inverted  as  they  were  for  the  case  of  (natural) 
thermal  emission.  This  latter  circumstance  is  due  primarily  to  the  absence  of  the  prefactor  (l-co0) 
in  the  definition  of  eq.  (29b). 


In  the  lower  plots  of  figure  8  the  abscissa  for  both  cases  is  the  slant  path  angle  cosine,  p,  and,  as 
before,  the  various  curves  in  each  of  the  two  subplots  refer  to  cloud  albedo.  Note  again  that,  as  in 
the  upper  plots,  the  ordering  of  the  curves  for  reflection  and  emission  are  the  same;  that  is,  both 
increase  with  increasing  albedo.  The  smooth  curves  in  figure  8  were  all  generated  using  the 
PILOT-EX  methodology  as  modified  to  treat  embedded  sources.  The  discrete  data  points  for  the 
reflective  case  were  obtained  directly  from  the  Van  de  Hulst  tables  using  eq.  (28)  and  those  for 
emission  were  obtained  directly  from  the  plots  of  figure  9.14  in  the  Van  de  Hulst  text. 


For  the  reflective  cases  (lower  left  plots),  it  is  significant  that  all  curves  exhibit  a  monotonic 
decrease  in  reflectivity  as  the  zenith  angle  cosine  increases.  For  example,  in  the  conservative 
case,  the  curve  begins  with  a  value  identically  equal  to  unity  at  u=0  (horizontal  propagation)  and 
eventually  decreases  to  a  value  of  around  0.610  at  u=l  (vertical  propagation).  This  general 
behavior  is  characteristic  of  all  curves  and  should  be  expected  since,  according  to  eq.  (26),  the 
slant  path  optical  thickness  (x/u)  approaches  infinity  as  u  approaches  zero  and  decreases  as  the 
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slant  angle  decreases.  As  the  albedo  decreases,  the  curves  all  decrease  in  magnitude  and  tend  to 
“flatten”  and  become  nearly  independent  of  angle  at  the  lowest  albedo. 


For  the  emissive  cases  (lower  right  plots),  the  curves  are  similar  to  the  reflective  case  in  general 
behavior  except  for  the  curious  maxima  that  occur  in  the  proximity  of  (1=0.2  [0=78°  ]  which  is 
most  pronounced  at  the  higher  albedos.  Again,  as  the  albedo  decreases  the  curves  tend  to 
decrease  in  magnitude  and  flatten  out  as  was  also  seen  for  reflective  case.  Although  it  may  not  be 
obvious  from  the  plots  it  does  turn  out  to  be  a  fact  that  the  embedded  source  results  can  be 
related  to  emissivity  through  a  scaling  factor  of  the  form  (l-co0).  We  have  made  use  of  this  fact  in 
table  2  where  we  again  summarize  results  in  numerical  form  for  both  the  reflective  and  emissive 
cases. 

Table  2.  Comparison  of  calculated  reflectivities  and  emissivites  for  various  slant  angles. 


REFLECTION 


R (t, 1) 

wo.  99 

wo=0 .80 

wo=0 .40 

wo=0 .20 

u 

Ref.  1 

Ref.  5 

Ref.  1 

Ref.  5 

Ref.  1 

Ref.  5 

Ref.  1 

Ref.  5 

.00010 

.97924 

.97960 

.66525 

.66973 

.26169 

.26639 

. 12012 

.12286 

.10000 

.97676 

.97675 

.63523 

.63518 

.23206 

.23272 

.10331 

.10329 

.30000 

. 94032 

.94032 

.59302 

.59300 

.20536 

.20534 

.08941 

.08940 

.50000 

. 84281 

. 84281 

.52782 

.52778 

.  18048 

.18046 

.07815 

. 07814 

.70000 

.74104 

.74103 

.46300 

.46298 

.15765 

.15763 

.06814 

.06813 

.  90000 

.65373 

.65373 

.40804 

.40803 

.13868 

.13867 

.  05989 

.  05988 

1.00000 

.61602 

.61602 

.38439 

.38437 

.13056 

.13055 

. 05637 

.05636 

EMISSION 


E(t,l) 

wo=0 .99 

wo=0 .80 

wo=0 .40 

wo=0 .20 

u 

Ref.  1 

Ref.  5 

Ref .  1 

Ref.  5 

Ref.  1 

Ref.  5 

Ref.  1 

Ref.  5 

.00010 

.02075 

.02040 

.33475 

.33027 

.73831 

.73361 

. 87988 

.87714 

.10000 

.02320 

.02320 

.36473 

.36477 

.76790 

.76723 

.89664 

.89666 

.30000 

.02400 

.02401 

.37130 

.37133 

.75896 

.75899 

.87492 

. 87493 

.50000 

.02185 

.02185 

.33685 

.33688 

.68418 

.  68420 

.78652 

.78652 

.70000 

.01931 

.  01932 

.29734 

.29737 

.60270 

.60272 

.69221 

.69222 

.90000 

. 01707 

. 01708 

.26276 

.26278 

.53213 

.53214 

.61092 

.61093 

1.00000 

.01610 

.01610 

.24773 

.24775 

.50156 

.  50157 

.57575 

.57576 

In  table  2,  the  upper  set  refers  to  reflectivity  and  lower  set  refers  to  emissivity  using  either  our 
methodology  and  eq.  28  or  the  angular  dependent  tables  of  Van  de  Hulst.  The  study  variables  in 
this  case  are  the  slant  angle  and  albedo  for  a  single  cloud  optical  thickness  of  unity.  Examination 
of  the  results  again  show  near  perfect  agreement  with  the  largest  errors  occurring  at  the  lowest 
angles  which  is  a  result  we  may  expect  due  to  the  very  high  slant  value  optical  thickness,  which 
complicates  the  multiple  scattering  solutions. 

Taken  together,  the  comparisons  here  for  both  reflection  and  emission  are  favorable  and  strongly 
support  the  validity  of  our  approach.  We  have  other  anecdotal  and  qualitative  evidence  from  the 
literature  on  both  cylinders  and  optically  thick  spheres  that  also  tend  to  support  our  conclusions, 
but  nothing  as  extensive  as  the  Van  de  Hulst  study. 
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6.  Bandpass  Functions  and  Apparent  Temperature 


It  has  been  assumed  throughout  that  the  bandpass  functions  introduced  in  [1]  take  care  of  all  the 
conversions  from  temperature  to  (thermal)  radiance  and  vice  versa.  That  is,  for  any  real, 
calibrated,  radiometric  measurement,  in  the  absence  of  any  emissive  or  obscuring  effects,  the 
sensor  signal  responds  in  accordance  with  the  following  expression: 


JUW,T)dX 

B»(T)  =  ^ - <— )  =  F;(T)^-  (30) 

|B(X,T)dX  71  n 

0 


where  X  is  the  wavelength  within  the  specified  window,  AX,  and  fAx(X)  represents  some  specified 
instrument  bandpass  “filter”  function.  Other  quantities  are  the  Stefan  Boltzman  constant,  o,  the 
(absolute)  real  thermodynamic  temperature  of  the  target  surface,  T,  and  the  Planck  Blackbody 
function,  B(X,T).  For  an  infrared  sensor,  eq.  (30)  produces  a  one-to-one  correspondence  between 
sensor  signal  and  surface  temperature,  provided  that  the  underlying  “target”  is  a  Blackbody  of 
uniform  temperature  and  there  are  no  intervening  emissive,  atmospheric,  or  obscurant  effects. 

To  illustrate,  in  figure  9,  we  show  plots  of  eq.  (30)  for  several  idealized  bandpass  functions  of 
practical  interest,  including  the  full  spectrum  case,  [i.e.,  fA>.(  X)  3  for  all  X],  for  which  the 
expression  reduces  to  the  familiar  Stefan-Boltzman  Equation  [i.e.,  B^  (T)  =  aT4/jr]. 
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—  full  spectrum 


—  full  spectrum 


Figure  9.  Bandpass  functions:  radiance,  left;  thermal  derivative,  right. 

In  making  the  plots  of  figure  9  we  have  modeled  the  filter  function,  fAx(X)  as  an  idealized  unit 
rectangular  step  function  over  the  specified  window  (AX),  as  indicated  on  the  figures.  The  plots 
on  the  left  are  the  bandpass  functions  and  those  on  the  right  are  the  corresponding  “thermal” 
derivatives  used  in  defining  the  “radiation  contrast”  that  will  be  explained  more  fully  in  a  later 
section.  In  all  cases,  the  mathematical  calculations  were  performed  using  the  semi-analytical 
methods  described  in  the  SPIE  (International  Society  for  Optical  Engineering)  Handbooks  [9] 
and  are  accurate  to  within  O.Olpercent  for  the  range  of  temperatures  shown.  For  ease  of 
explanation  we  have  assumed  that  the  idealized  instrument  has  been  calibrated  against  a 
Blackbody  target  so  that  any  instrumental  artifacts  are  either  compensated  or  otherwise  removed 
[10].  Thus,  with  an  instrument  so  calibrated,  it  is  possible,  in  principle,  to  convert  the  measured 
signal  to  the  target  temperature  using  the  inverse  of  eq.  (30),  or  graphically  using  the  plots. 
Conversely,  if  we  know  the  target  temperature,  we  can  immediately  determine  the  corresponding 
radiance  using  either  the  graphs  or  eq.  (30)  in  the  forward  mode.  Thus  we  can  think  of  the  plots 
as  a  “calibration”  function  for  converting  radiometric  measurements  to  temperature  and  vice 
versa  as  long  as  the  “targets”  are  blackbodies  and  there  are  no  intervening  atmospheric  or  aerosol 
effects. 


To  include  the  effect  of  a  “hard  target”  emissivity,  the  conventional  approach  is  to  define  an 
“apparent”  temperature  [11]  such  that  the  total  energy  in  the  bandpass  is  the  same  as  that  in  the 
reduced  signal.  That  is,  we  form  the  following  “energy  balance,”  accounting  for  both  direct 
emissions  from  the  target  and  any  reflections  from  the  ambient  surroundings,  as 
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(31) 


JfM(\)B<3.,T; )  =  e,„  JfM  (X)B(X,T„  )dX 

0  o 


+  V  Jfa(X)B(X,T;„)dX 


where  T*tgt  is  defined  as  the  “apparent”  target  temperature  and  Ttgt  is  the  corresponding  real  (i.e., 
thermodynamic)  target  temperature.  Other  terms  in  eq.  (31)  are  the  target  emissivity,  etgt,  the 
target  reflectivity,  rtgt,  and  the  apparent  ambient  temperature,  T  ext>  which  determines  the 
magnitude  of  incoming  irradiance  from  the  ambient  (sky  and  surface)  surroundings.  Note  that 
here  and  throughout  we  adopt  a  notation  wherein  the  superscript  (*)  denotes  apparent 
temperature. 

Ignoring  reflections  for  the  moment,  it  is  clear  that  the  effect  of  emissivity  is  to  cause  a  uniform 
signal  reduction  across  the  spectrum  and  is  quite  different  than  that  of  changing  the  temperature 
which  actually  causes  a  non-uniform  change  in  the  underlying  spectrum  as  well  as  a  reduction  in 
the  overall  signal  magnitude.  The  situation  is  illustrated  in  figure  10  where  we  have  plotted  a 
series  of  spectra  corresponding  to  various  apparent  temperatures  on  the  left  and  the 
corresponding  energy  equivalent  spectra  of  the  emissive  counterparts  on  the  right.  For 
convenience  we  use  the  “full  spectrum”  relationship  (i.e.,  Tapp  =  e1/4  Treai). 


A  subtle  but  significant  difference  in  the  two  sets  of  curves  is  the  shift  of  the  maximum  to  longer 
wavelengths  in  the  left  set  but  a  simple  multiplicative  effect  for  each  curve  on  the  right  set.  This 
inconsistency  is  a  consequence  of  the  definition  of  apparent  temperature  that  can  cause  problems 
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in  applications  where  it  is  critical  to  determine  the  actual  value  of  the  (real)  temperature.  For  our 
application  however,  we  are  ultimately  concerned  only  with  measured  radiant  energies  for  which 
this  inconsistency  is  not  a  major  concern. 

The  question  now  arises  as  to  how  to  treat  the  received  signal  in  the  presence  of  an  obscuring 
aerosol  which  we  assume  can  either  reduce  or  enhance  the  received  signal  due  to  the  actions  of 
extinction,  emission,  and  ambient  in-scatter  plus  the  possibility  that  the  cloud,  background,  and 
target  temperatures  may  all  be  different.  Ignoring  for  the  moment,  the  effects  of  emission  and 
scattering,  it  is  clear  that  the  effect  of  direct  (Beers’s  Law)  signal  attenuation  can  be  treated  by 
extending  the  usual  definition  to  include  a  further  signal  reduction  in  the  “energy  balance” 
relationship  as  follows: 


Jftt(X)B[X,T;(i)]  =  t,„W{Jf4l(X)8wB(A,Tw)ca  +  Jf^Wr^BPXjdX}  (32) 


where  taer  (x)  [=e'T]  is  the  cloud  transmissivity  and,  T*tgt  (x)  is  the  target  apparent  temperature  as 
sensed  through  an  obscurant  of  optical  thickness  (x)  and,  for  the  moment,  accounting  only  for  the 
direct  attenuation.  Thus,  the  treatment  for  transmittance  is  similar  to  that  for  “hard  target” 
emissivity  and  the  combined  effect  of  both  is  multiplicative.  Note,  too,  that  this  “apparent”  target 
temperature  is  the  radiative  temperature  that  would  actually  be  sensed  by  an  instrument 
calibrated  against  a  Blackbody. 

The  inclusion  of  the  effects  of  the  aerosol  emissivity  and  reflectivity  is  treated  independently,  but 
in  a  similar  maimer.  That  is,  accounting  for  the  total  diffuse  (path)  radiance  from  the  cloud,  we 
form  the  appropriate  energy  balance  to  define  the  apparent  cloud  temperature,  T  cid(x),  analogous 
to  the  apparent  target  temperature,  as  follows: 


JfMWB[X,T;u(T)]  =  e„(T)JfM(MB(>..Tck,)dX  +  r>tr(T)|fal(>.)B(>.,T;,)dX  (33) 

0  0  0 


where  eaer(x)  and  raer(x)  are,  respectively,  the  aerosol  emissivity  and  reflectivity  as  defined  in 
previous  sections,  Tcid  is  the  real  (thermodynamic)  temperature  of  the  cloud,  and  Text  remains  as 
defined  in  eq.  (31).  Note  that  the  above  expression  defines  the  apparent  cloud  temperature  only, 
not  including  the  “target”  contribution. 

We  now  define  the  composite  apparent  temperature  as  affected  by  all  processes  by  adding 
contributions  from  both  the  (attenuated)  target  signal  (eq.  32)  and  cloud  radiance  (eq.  33)  and 
thus  obtain,  for  the  total  energy  balance,  the  following: 


Jf41(X)Bp.,T;„(T)]=  Jf4l(X)B[X,T;(T)]d3.+  Jf^POBtWlWldX 


(34) 


where  we  explicitly  assert  that  the  apparent  temperature  in  this  case  represents  the  actual 
“measured”  radiometric  temperature,  Tmea(x),  as  would  be  sensed  by  a  real  instrument.  This 
assertion,  which  may  appear  at  first  glance  as  a  leap  of  faith,  is  clarified  by  realizing  that  the 
whole  derivation  is  based  on  a  careful  treatment  of  the  energy  balance  at  each  step  using  what  we 
assume  to  be  “measured”  data. 

Up  to  this  point,  we  summarize  by  expanding  and  rewriting  eq.  (34)  more  conveniently  in  terms 
of  the  bandpass  functions  of  eq.  (30)  as  follows: 

Blx  [T»  (*)]  =  t  aer  (x)B*  <T; )  +  e  aer  (x)B*  (Tcld )  +  raer  (x)B*  (Tc*xt )  (35) 


where  T*tgt  is  the  target  apparent  temperature  determined  via  eq.  (31).  At  this  point  it  should  be 
made  clear  that,  for  a  given  instrument  specification,  all  terms  on  the  right  side  of  eq.  (35)  can  be 
obtained  directly  from  calibration  curves  such  as  those  of  figure  10  (or  equivalent  from 
numerical  look  up  tables)  using  the  known  temperatures  and  the  calculated  aerosol  properties. 
The  sum  of  all  terms  then  gives  the  full  signal  radiance  from  which  the  apparent,  or  “measured,” 
temperature  can  be  obtained  by  inverting  the  process  with  the  same  calibration  curve. 
Alternatively  we  can  express  this  intuitive  process  formally  as, 

t„  =B2;’>(Rn,„) 

where  (36) 

taer  (x)B*  (T; )  +  {Eaer  (x)B*  (Tcld )  +  raer  COB*  (T*xt)} 


where  the  superscript  (-1)  represents  the  inverse  of  the  function  and  all  other  quantities  have 
been  identified  in  the  previous  paragraphs.  As  usual,  the  first  term  on  the  right  side  of  eq.  (36)  is 
the  direct  contribution  and  the  second  (bracketed)  term  is  the  diffuse  contribution  (also  called  the 
path  radiance  in  applications). 

We  now  work  an  example  using  the  hypothetical  data  of  table  3  assuming,  for  convenience,  a 
“full  spectrum”  bandpass  and  a  cloud  total  optical  thickness  of  unity  (i.e.,  taer=e'!  =  0.368). 


Table  3.  Hypothetical  data  for  bandpass  example. 


EXAMPLE  r=1.0 

units  =  watt/m2  Sr'1 

ta„=0.368 

Tcld=40°  C 

Bcid=173 

W*BcH=63.9 

£^=0.421 

t\x,=io°c 

B„,=116 

£aer*Bext=48.8 

raer=0.211 

T*i«=20°C 

B,m=133 

r«r*Btrt=28.1 

I^mea 

=140.8 

In  table  3,  the  first  two  columns  of  data  are  the  raw  input  that  we  have  arbitrarily  “made  up”  for 
the  example  and  the  data  in  the  remaining  two  columns  are  calculated  using  the  methods 
described  above.  The  first  step  in  the  example  is  to  calculate  the  fractional  radiance,  B*ax  (T),  for 
each  of  the  three  components.  The  computations  can  be  done  either  graphically  from  the  upper 
leftmost  plot  of  figure  9  or,  since  we  are  addressing  the  full  spectrum  bandpass,  from  the  Stefan- 
Boltzman  analytical  expression.  In  either  case  the  results  are  as  shown  in  the  third  column  of 
table  3.  The  data  in  fourth  column  represent  the  products  of  column  1  and  column  3  which  give 
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each  of  the  three  terms  on  the  right  side  of  eq.  (36).  The  resultant  total  measured  radiance  is 
found  by  summing  the  final  column,  which  yields  a  value  of  140.8  (watt/m2  Sr'1 )  from  which  we 
derive  an  apparent  temperature  of  24.1  °C  by  the  inverse  procedure.  Note  from  the  example  that 
if  we  had  ignored  the  path  radiance  contributions,  the  resultant  radiance  of  63.9  watt/m  Sr' 
would  correspond  to  an  apparent  temperature  of -29. 2  °C,  producing  an  error  of -73.3  °C,  which 
is  significant  in  most  applications. 


7.  Radiation  Contrast  (Thermal  Derivative) 


It  turns  out  that  scanning  systems  such  as  a  FLIR  (forward-looking  infrared)  imager  are  not 
designed  to  rely  upon  the  magnitude  of  a  sensed  signal,  per  se,  but  instead  utilize  the  concept  of 
“ac  coupling”  which  is  based  upon  sensing  incoming  signal  differences,  or  contrast.  In  the 
analysis  of  these  type  systems  it  is  usual  to  address  the  case  wherein  the  sensed  temperature 
differences,  AT,  are  small  enough  so  that  we  can  use  the  following  linear  expansion: 

dR 

I(T  +  AT)  =  I(T)  +  ^-AT  (37) 

dT 


where  I(T)  and  I(T+AT)  represents  two  measured  signals  received  from  adjacent  scene  “pixels” 
and  R[=R(AX,,T)]  will  later  be  identified  with  the  band  weighted  blackbody  function  of  eq.  (30). 
Upon  rearranging  eq.  (37),  we  come  about  the  definition  for  radiation  contrast  [9],  Cr(T),  or  what 
actually  amounts  to  the  thermal  derivative,  that  is  [11], 


Cr(T) 


I(T  +  AT)  -  I(T)  _  AI  _  dR 
AT  ~  AT  ~  dT 


(38) 


where  we  have  taken  care  to  point  out  and  make  explicit  that  the  thermal  derivative  is  to  be 
evaluated  at  the  measured  temperature,  Tmea  (i.e.,  the  apparent  “temperature”  that  the  receiver 
actually  senses),  and  must  include  the  full  effect  of  any  intervening  aerosols,  including  path 
radiance.  For  our  application,  the  thermal  derivative  can  be  evaluated  directly  from  eq.  (36);  that 
is,  generalizing  somewhat  we  have 


taerCO 


dR, 


dT 


dR2 

+  Eaer(t)— 


dT 


,  ^R3 
+  raer(T> 


dT 


(39) 


where,  for  our  application,  we  will  eventually  replace  the  R  functions  with  the  fractional 
blackbody  functions  of  eq.  (30).  Our  point  in  writing  eq.  (39)  as  such  is  to  make  it  clear  that  if 
the  functional  forms  of  Ri,  R2,  and  R3  are  the  same,  then  the  three  derivatives  are  all  the  same 
since  they  are  all  evaluated  at  the  same  point,  Tmea.  Thus,  when  actually  applied  to  eq.  (36)  we 
have 
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MtacrOO  +  SaerCO  +  raerCO) 


(40) 


dR 

dT 


dBlx(T) 


dT 


dBlx(T) 


dT 


where  in  writing  the  last  expression  we  have  used  the  continuity  condition  [e(x)+r(x)+t(x)=T]  of 
eq.  (25)  and,  as  before,  B*(T)  is  the  bandpass  function  appropriate  for  the  particular  system 
under  consideration.  Some  plots  of  the  radiation  contrast  for  various  bandpasses  of  interest  were 
presented  earlier  in  the  rightmost  plots  of  figure  9.  It  is  worth  emphasizing  that  eq.  (38)  applies 
to  the  situation  where  the  contrast  is  small  and  should  not  be  used  for  situations  where  the 
contrast  is  large  (in  which  case  the  target  is  usually  obvious  anyway). 

The  impact  of  the  radiation  contrast  on  the  actual  “sensed”  temperature  difference  can  be 
determined  by  combining  the  definition  of  eq.  (38)  with  eq.  (40)  to  obtain 


dT 


(AT)mea  =CR(Tmea  )(AT)n 


(41) 


where  (AI)mea  represents  the  measured  ac  signal  difference,  and  (AT)mea  is  the  corresponding 
equivalent  temperature  difference.  It  is  also  useful  to  apply  a  form  of  eq.  (41)  to  the  system  noise 
(radiance)  signal,  ARsys,  to  produce 


ARsys=CR(Tn)ATn  (42) 

where  ATn  is  called  the  “noise  equivalent”  temperature  differential  commonly  used  in  hardware 
systems  analysis. 

We  now  see  an  immediate  problem  for  the  analyst  because  the  exact  measured  temperature,  Tmea, 
is  not  necessarily  known  for  these  types  of  system  since  only  differences  are  considered. 
Nevertheless,  the  concept  of  sensing  contrast  directly  is  an  attractive  alternative  in  applications 
because  there  is  at  least  a  partial  mitigation  of  the  (usually  unfavorable)  effects  of  path  radiance 
which  comes  about  by  differencing  (see  section  7).  It  is  interesting  to  note,  however,  that  the 
effect  of  path  radiance  on  radiation  contrast  can  be  to  enhance  the  sensitivity.  This  can  be  seen 
from  the  rightmost  plots  of  figure  9  where  it  is  clear  that  the  radiation  contrast,  and  thus  the 
system  sensitivity,  generally  increases  as  the  sensed  temperature  increases.  However,  this  is  not 
the  complete  story  since  we  have  yet  to  examine  the  effect  of  aerosol  induced  noise,  which  we 
will  consider  in  the  next  section. 
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8.  Effect  of  Aerosol-Induced  Noise 


As  noted  in  the  previous  section,  in  modeling  the  performance  of  an  ac  coupled  system  such  as 
the  FLIR,  the  received  radiance  can  be  modeled  as  the  difference  between  two  successive 
samples  from  adjacent  pixels,  and  one  of  the  relevant  parameters  affecting  system  performance  is 
the  signal-to-noise  ratio,  which  we  will  derive  in  this  section.  In  this  case  the  incident  “signal”  in 
the  presence  of  an  intervening  aerosol  is  modeled  as 

AI(Xj  ,x2 )  =  [R(x, ,  T,  mea )  -  R(t2  ,  T2>mea )]  (43) 


which  represents  the  difference  in  the  measured  radiance  between  two  adjacent  “pixels”  labeled 
by  subscripts  “1”  and  “2.”  Applying  eq.  (36)  twice  for  the  two  pixels  in  question  immediately 
yields 


AICej  5  x2 )  =  [e_t‘  (T*tgt )  -  e"*2  (T2  ,gt )]  +  [R(Xj ,  T,>dif )  -  R(t2  ,  T2jdif )] 


(44) 


where  we  have  arranged  terms  such  that  the  first  bracketed  expression  is  the  direct  contribution 
to  the  differential  signal  (i.e.,  ARdir)  and  the  second  represents  the  aerosol  diffuse  contribution 
(i.e.,  ARdif),  also  called  the  path  radiance. 

We  next  follow  the  usual  development  by  invoking  the  following  approximations: 


T,  ST2 


T  —  T  ~  T 

xl,dif  —  x2,Tdtf  Adif 


(45) 


both  of  which  may  usually  be  justified  because  the  two  pixels  are  adjacent  so  that  the  optical 
paths  through  the  atmosphere  are  nearly  the  same.  Applying  the  above  approximation  to  eq.  (44) 
yields 


AI(t,  ,  t2  )  =  e~T  [B^  (Tj*tgt )  -  (T2  tgt )]  +  [R(Xj ,  T,  dif )  -  R(t2  ,  T2  dif )] 

=  [e'TB*AX(Ti;tgt)-B*AX(T2*|tgt)]  +  [ARaer(x,Tdif )] 


(46) 


where  the  second  bracketed  expression,  ARaer(x,  Tdif ),  represents  the  path  radiance  (difference) 
contribution  which,  owing  to  the  approximations  of  eq.  (45),  is  assumed  to  be  small  but  not 
necessarily  zero  and  not  necessarily  small  with  respect  to  the  direct  term.  It  is  natural  to  take  the 
first  bracketed  term  in  eq.  (46)  to  be  the  true  difference  signal  (i.e.,  derived  from  that  part  of  the 
received  radiance  that  actually  originated  from  the  target)  and  the  remainder  is  taken  as  noise. 
Thus,  the  signal  to  noise  ratio  becomes 
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SNR(x)  = 


e'PLcr^-Bl.cr^)] 

ARaer(T,Tdif)  +  ARsys 


(47) 


where  we  have  also  added  a  term,  ARsys,  to  include  any  noise  contribution  from  the  system 
hardware  as  is  common  in  “clear  air”  studies  more  often  addressed  in  the  literature  [11].  It  is 
often  tempting  to  assume  that  the  path  radiance  contributions  either  cancel  each  other,  or  at  least 
the  difference  is  small  compared  to  the  direct  term  in  which  case  the  aerosol  noise  term  can  be 
neglected.  At  certain  extremes  this  is  a  valid  assumption  and  a  case  that  we  will  address  later; 
however,  in  general,  it  is  more  reasonable  to  assume  the  aerosol  contribution  to  be  proportional 
to  the  (average)  value  of  the  magnitude  of  the  path  radiance.  Thus  for  later  reference  we  write 


AR  (T  T  OC  (^•(Tl>^l.llif)  +  ^(T2>rr2,dif)  _  c  r>  '"T  \ 

^■aer(T’  -*■  dif  )  K  „  ~dwR(X,Tdif) 


(48) 


where  8W  is  an  ad  hoc  “weather”  dependent  parameter  which  we  conjecture  can  range  from  a 
value  of  zero  for  low  turbulence  conditions  to  a  value  on  the  order  of  unity  for  high  turbulence 
conditions.  Clearly,  the  path  radiance  in  eq.(48)  can  be  written  in  terms  of  the  aerosol  cloud 
properties  as 


R(tX  )  =  eaer(x)B^(Tcld)  +  U  t)B;xCO  (49) 

where  Tcid  is  the  cloud  temperature,  and  T*ext  is  the  ambient  radiative  temperature,  as  defined  in 
connection  with  eq.  (36)  of  the  previous  section. 

Continuing  with  the  usual  development,  we  next  invoke  the  “small  signal”  approximation  which 
amounts  to  assuming  that  the  underlying  direct  signal  difference  is  small  enough  such  that 
eq.  (38)  of  the  previous  section  holds.  In  this  case  we  can  replace  the  radiance  differentials  in  the 
numerator  of  eq.  (47)  with  their  corresponding  apparent  temperature  differentials  as  follows: 

Kx  CC )  -  BL (T2*,lgt )  *  cr (Tmta  )(T*tgt  -  T2*tgt )  (50) 


where  Cr  (Tmea )  is  the  thermal  derivative  defined  by  eq.  (41)  in  the  previous  section  and  is 
evaluated  at  the  measured  apparent  temperature  which  we  will  ultimately  take  as  the  average  of 
the  two  target  pixels. 

Making  the  obvious  substitutions  in  eq.  (47),  we  now  have 


SNR(x)  = 


e~TCr(T0)(AT)0 
5wR(t,Td*if)+  Cr(Tn)ATn 


(51) 


where  we  have  also  used  the  thermal  derivative,  Cr(Tn),  to  convert  the  system  noise,  Rsys,  from  a 
radiance-based  parameter  to  a  noise  equivalent  temperature  difference,  ATn,  as  it  is  commonly 
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done  in  systems  analysis  [cf.,  eq.  (42)].  We  have  also  substituted  T0  for  the  average  apparent 
temperature  [i.e.,  T0=(Ti>tgt  +T2>tgt)/2]  and  AT0  as  the  corresponding  temperature  differential  (i.e., 
AT0=Ti,tgt  -T2,tgt). 


Finally,  collecting  all  terms  and  making  all  of  the  suggested  substitutions,  we  have 


SNR(x)  =  { 


8wR(T,Td‘if)  +  Cr(Tn)ATn 


Cr(T0)(AT)0 

B^cr.) 


(52) 


where,  for  convenience,  we  have  multiplied  both  numerator  and  denominator  by  the  average 
target  radiance,  B*ax  (T0).  Written  as  such  the  second  bracketed  term  in  eq.  (52)  represents  the 
surface  level  (i.e.,  unobscured)  radiance  based  contrast  at  the  target  surface  (i.e.,  ARo/Ro)  and  the 
first  bracketed  term  represents  a  generalized  version  of  what  is  sometimes  called  the  “contrast 
transmission,”  the  reciprocal  of  which,  for  historical  reasons,  is  referred  to  in  the  older  literature 
as  the  “sky-to-ground”  ratio  (although  the  terminology  is  somewhat  inappropriate  here). 

In  the  limit  where  the  aerosol  noise  can  indeed  be  neglected  either  because  conditions  are  such 
that  the  parameter,  8W,  is  small  or  because  the  system  noise,  ARsys,  is  relatively  large,  eq.  (52) 
reduces  to  the  commonly  accepted  “clean  air”  form  given  by 

SNR(T)4K^=e-’-^  (53) 

AT„ 


where  we  have  assumed  that  the  temperatures  for  the  signal  and  system  noise  are  sufficiently 
close  that  the  thermal  derivatives  are  the  same  (i.e.,  T0~Tn),  an  assumption  which  requires  some 
care  and  attention  in  real  applications. 

Another  example  for  which  we  have  some  interest  is  the  case  of  thermal  equilibrium  between  the 
cloud  and  the  ambient  atmosphere  (i.e.,  Tcid=Text )  in  which  case,  neglecting  the  system  noise 
term,  the  contrast  transmission,  Tcon,  becomes 


X  _  r  e  i  r  ®AxCf0)  ^ 

T-n  (T)  -  t "I 


1-e" 


B  a(TMy 


(54) 


where  we  have  also  used  the  general  condition  for  isothermal  clouds  [e(T)+r(x)=l-e‘T  ]  from 
eq.  (25).  For  the  full  equilibrium  case  (i.e.,  Ttgt=Tcid=T0),  eq.  (54)  reduces  to  the  well-known 
classical  result  [(l-e~T)/e'T],  which  is  often  used  for  plane  layers.  Thus  it  appears  that  the 
development  here  ties  together  several  often  used  approximations  and  offers  a  way  to  bridge  the 
various  extremes.  The  problem  at  the  moment  is  that  the  form  of  the  ad  hoc  weather  dependent 
parameter,  8W,  is  not  known  and  does  not,  at  this  period  in  time,  appear  to  be  tractable 
theoretically.  We  do  have  an  ongoing  experimental  effort  based  on  (near)  simultaneous  field 
measurements  of  path  radiance  and  transmittance  to  shed  some  light  on  this  issue  [13]. 


29 


9.  Temperature  Stratification 


In  this  section  we  extend  the  study  to  include  temperature  stratification  in  the  emissive  cloud 
formulation.  The  extension  is  reasonably  straightforward  except  that  the  definition  of  the  cloud 
emissivity  is  no  longer  consistent  with  basic  physical  principles  as  it  was  for  the  isothermal  case. 
However,  this  is  only  a  minor  problem  that  we  circumvent  by  defining  a  “practical”  emissivity 
normalized  to  the  maximum  which,  if  used  consistently,  causes  no  real  problem  in  applications 
of  our  type.  Also,  for  the  non-isothermal  case,  the  emissive  source  functions  no  longer  mirror  the 
inverse  of  the  reflective  source  functions  as  they  did  for  the  isothermal  case  and  therefore  the 
general  expression  relating  emissivity  to  reflectivity  [eq.  (25)]  is  no  longer  valid  and  other 
fundamental  notions  need  to  be  modified.  However,  no  real  conceptual  difficulties  exist  in 
actually  applying  the  method,  but  the  treatment  is  somewhat  less  elegant  and  the  results  a  bit 
more  complex. 

As  before,  the  underlying  equation  is  the  thermal  version  of  the  radiative  transfer  equation 
which,  for  reasons  cited  above,  we  normalize,  in  first  order,  as  follows: 


(1-cQq) 

Blx(T«) 


(55) 


where  t0  and  to0  are,  respectively,  the  cloud  total  optical  thickness  and  albedo,  x  is  the  optical 
thickness  referenced  to  zero  at  the  edge  of  the  cloud  (see  sec.  2)  and  T^x  is  the  maximum  cloud 
temperature  which,  for  our  hotter-than-ambient  Gaussian  clouds,  will  come  about  at  the  cloud 
center. 

The  cloud  temperature,  Tcid  (r),  is  the  aerosol  real  (thermodynamic)  temperature  at  a  radial 
distance,  r,  from  the  cloud  center  and  is  modeled  in  our  standard  Gaussian  cloud  as 

Tcid (r)  =  Tamb  +  Pmix (r)0™  - Tamb )  (56) 

where  Tamb  is  the  ambient  temperature  outside  the  cloud  and  pmjx(r)  is  a  “mixing  factor”  taken  as 
proportional  to  the  aerosol  concentration  normalize  to  unity  at  the  cloud  center  [i.e.,  pmjX  (r)= 
C(r)/C0]  and  where  the  aerosol  concentration,  C(r),  is  modeled  in  accordance  with  the  procedures 
described  in  section  2.  It  is  clear  from  inspection  that  eq.  (56)  yields  a  maximum  value,  Tmax,  at 
the  cloud  center  and  drops  to  near  ambient,  Tamb,  at  the  cloud  edge.  It  is  furthermore  clear  that 
the  overall  radial  symmetry  is  preserved  so  that  the  multiple  scattering  methodology  described  in 
our  previous  report  remains  valid;  thus,  we  can  calculate  the  multiple  scattering  source  function 
iteratively  using  eq.  (55)  as  the  initial  condition.  Some  examples  of  thermally  stratified  emissive 
source  functions  are  shown  in  figure  11  for  a  variety  of  cloud-to-ambient  temperature 
differentials. 
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In  figure  1 1,  the  first  row  of  subfigures  represents  the  isothermal  case  which  is  identical  to  the 
corresponding  plots  of  figure  7  but  which  we  repeat  here  for  reference.  The  remaining  three  rows 
correspond  to  the  various  temperature  differentials  [AT=Tmax-Tamb]  as  indicated  to  the  right  of 
each  row.  Each  of  the  four  columns  correspond  to  a  different  cloud  optical  thickness  ranging 
from  xo=0.50  to  4.0  as  indicated  at  the  top  of  each  column.  Each  subplot,  of  which  there  are  16 
total,  includes  1 1  curves  (9  discemable)  corresponding  to  different  values  of  albedo  increasing  in 
equal  increments  from  ©o=0  (extreme  top)  to  co0=l  (extreme  bottom)  following  the  usual 
monotonic  scheme  used  throughout. 


In  studying  the  rows  and  columns  of  figure  11,  two  general  trends  are  immediately  apparent  as 
manifested  by  the  increased  curvature,  or  “narrowing,”  of  the  curves  as  either  the  optical 
thickness  increases  (along  each  row)  or  as  the  temperature  differential  increases  (down  each 
column).  There  is  also  a  small  but  noticeable  increase  in  the  magnitude  of  the  source  function 
with  increasing  optical  thickness,  which  seems  to  be  most  marked  at  the  higher  albedos 
represented  by  the  lower-lying  curves  in  each  subplot.  The  effect  of  increased  temperature  has 
little  effect  on  the  magnitude  of  the  source  functions  except  at  the  cloud  edges  where  the 
magnitude  clearly  decreases  with  increasing  temperature  differential,  a  behavior  due  primarily 
the  normalization  convention.  The  results  support  the  assumption  of  a  strong  build  up  of  internal 
radiance  in  the  center  as  the  temperature  differential  is  increased.  From  the  trends  shown,  we 
anticipate  that  the  curves  continue  to  narrow,  as  the  temperature  differential  increases,  and 
eventually  collapse  to  near  delta  functions  at  the  extreme,  especially  at  the  higher  values  of 
optical  thickness. 
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The  general  behavior  of  the  source  functions,  most  noticeably  the  build  up  of  source  radiance  in 
the  cloud  interior,  suggests  that  the  effect  will,  in  turn,  correspond  to  a  lowering  of  the 
emissivity.  This  is  intuitively  plausible  because  any  source  at  the  cloud  center  will  necessarily 
experience  some  “self  extinction”  in  reaching  the  cloud  exterior.  Thus,  we  would  expect  the 
calculated  emissivities  to  become  lower  and  lower  as  the  temperature  differential  increases  and 
become  more  and  more  pronounced  as  the  optical  thickness  increases.  Some  of  our  preliminary 
calculations  have  borne  this  out;  in  particular  we  have  found  emissivity  ratios 
[£(AT=x)/e(AT=0)]  to  be  on  the  order  of  52  percent  for  an  optical  thickness  of  x=\  and  on  the 
order  of  40  percent  for  an  optical  thickness  of  t=4,  and  in  both  cases  an  albedo  of  ©o=0.50. 
Unfortunately,  we  have  not  yet  found  a  way  to  present  the  results  in  a  systematic  way,  but  this  is 
a  viable  goal  for  further  work.  We  have  also  looked  at  some  data  on  emissive  sources  and  do 
have  evidence  that  self  extinction  occurs  in  real  smokes  [15]. 


10.  Summary 


The  details  given  here  along  with  those  also  given  in  our  companion  paper  [1]  offer  a  complete 
solution  to  the  radiative  transfer  equation  for  finite  aerosol  clouds  of  radial  symmetry.  The 
method  addresses  both  aerosol  thermal  emission  as  well  as  reflection,  including  all  orders  of 
multiple  scattering.  The  method  represents  exact  solutions  that  compare  favorably  with  other 
studies  when  applied  to  plane  layers.  We  have  also  demonstrated  the  significance  of  the  results  in 
systems  evaluation  through  the  effect  on  the  thermal  derivative,  or  “radiation  contrast,”  used  in 
evaluating  infrared  scanning  systems  such  as  the  FLIR.  The  methodology,  per  se,  is  applicable 
only  to  isothermal  clouds  and  isotropic  scattering  although  these  are  not  necessarily  fundamental 
limitations  on  the  method.  In  another  preliminary  study  we  have  extended  the  method  to  include 
non-isothermal  clouds  and  other  practical  effects  using  various  approximations  [13,14],  Other 
efforts  are  underway  to  evaluate  the  results  against  field  measurements  [16]. 
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